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ABSTRACT 

i-S^ ■ The outer giant planets, Uranus and Neptune, pose a challenge to theories of planet 

formation. They exist in a region of the Solar System where long dynamical timescales 

2 ■ and a low primordial density of material would have conspired to make the formation of 
c/3 . such large bodies (~ 15 and 17 times as massive as the Earth, respectively) very difficult. 
^ ' Previously, we proposed a model which addresses this problem: Instead of forming in 

the trans-Saturnian region, Uranus and Neptune underwent most of their growth among 
^ ' proto-Jupiter and -Saturn, were scattered outward when Jupiter acquired its massive 

■ gas envelope, and subsequently evolved toward their present orbits. We present the 

results of additional numerical simulations, which further demonstrate that the model 
readily produces analogues to our Solar System for a wide range of initial conditions. 
We also find that this mechanism may partly account for the high orbital inclinations 
observed in the Kuiper belt. 

Subject headings: celestial mechanics — solar system: formation — Kuiper belt — planets 
and satellites: formation 



1. Introduction 



The growth of Uranus and Neptune in the outer Solar System is not readily accounted for 
by conventional models of planet formation. A low primordial density of planetesimals and weak 
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solar gravity would have made the process of accretion slow and inefficient. In direct N-body 
simulations of accretion among (approximately) Earth-mass bodies beyond 10 AU, performed with 
three different computer codes, little accretion is found to take place over timescales of 10^ years, and 
by extrapolation, over the age of the Solar System (Levison & Stewart 2001). Earlier simulations 
by Brunini k. Fernandez (1999) showed accretion of the ice giants in several xlO^ years with the 
same initial conditions, but later simulations, performed with an improved integrator, require that 
bodies be enhanced in radius by at least a factor of ten, relative to bodies having the density of 
Uranus and Neptune, in order to recover the previous result (Brunini 2000). Therefore, Uranus 
and Neptune are unlikely to have formed from a late stage of mergers among large protoplancts, 
analogous to the putative final phase of planet formation in the terrestrial zone (eg. Wetherill 1996, 
Chambers & Wetherill 1998). The oligarchic growth model, in which the principal growth mode 
is accretion of small planetesimals by a protoplanet, also produces timescales which are too long 
(Kokubo and Ida 2000, Thommes 2000), unless the feeding zones of Uranus and Neptune can be 
replenished quickly enough with low random velocity planetesimals from elsewhere in the nebula 
(Bryden, Lin and Ida 2000). 

Thommes, Duncan &: Levison (1999), hereafter TDL99, develop an alternative model to in-situ 
formation for the origin of Uranus and Neptune. Beginning with four or more planetary embryos 
of 10-15 in the Jupiter-Saturn region, they explore through N-body simulation the evolution of 
the system after one of these bodies (and in one case, two at the same time) accretes a massive gas 
envelope to become a gas giant. They find that the remaining giant protoplanets are predominantly 
scattered outward. Dynamical friction with the planetesimal disk subsequently recircularizes their 
orbits, which leads, in about half the simulations performed, to a configuration quite similar to 
the present outer Solar System, with the scattered giant protoplanets taking the roles of Uranus, 
Neptune and Saturn. These results suggest that Uranus and Neptune are actually potential gas 
giant cores which formed in the same region as Jupiter and Saturn, but lost the race to reach 
runaway gas accretion. 

Here, wc explore this model in more detail, and perform further simulations. Sections 2 and 3 
motivate our choice of initial conditions for the simulations. Section 4 discusses the mechanism for 
transporting a proto-Uranus and -Neptune to the outer Solar System. N-body simulation results 
are presented in Section 5. The effect on the asteroid and Kuiper belts is discussed in Section 6. 
We summarize and discuss our findings in Section 7. 

2. Available material in the Jupiter-Saturn region 

Hayashi (1981) estimated the minimum primordial surface density of solids in the outer Solar 
System to be 



<ymin{a) = 2.1{a/bAU)-^'^ g/crr?. 



(1) 
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The requirement that Jupiter and Saturn's cores be massive enough to have initiated runaway 
gas accretion suggests that they are ~ 10 in mass (Mizuno et al 1978, Pollack et al 1996). 
Interior models are consistent with such a core mass, but also allow a coreless Jupiter (Guillot 
1999). We assume in this work that both Jupiter and Saturn began as ~ 10 bodies; putting 
gas giant cores and ice giant cores on the same footing is necessary in the "strong" version of our 
model, though it is not essential to the basic mechanism; we discuss variations on the model in 
Section 7. 

A surface density of 2.7 g/cvo? was likely too low to form a ~ 10 body at 5 AU before the 
gas was removed from the protoplanetary disk. Lissauer (1987) finds that a surface density of 15-30 
g/cm^ is needed to allow formation of Jupiter's core on a timescale of 5 x 10^ — 10^ years, while 
the model of Pollack et al (1996), which includes concurrent accretion of solids and gas, produces 
Jupiter in less than 10^ years with 10 g/cm^. The formation of giant planet cores may have been 
triggered at least in part by the enhancement in the solids surface density beyond the "snow line" , 
where water goes from being a gaseous to a solid constituent of the protoplanetary disk. In fact, 
outward diffusion and subsequent freezing of water vapor from the inner Solar System may have 
resulted in a large local density enhancement around 5 AU, perhaps yielding a surface density even 
higher than 30 g/cm^ (Stevenson &; Lunine 1988). Here we assume a power law surface density, 

(T{a)=(To{a/bAUy. (2) 

The above discussion suggests a ~ 10 to 30 g/cm^ as a plausible surface density at 5 AU. Allowing, 
also, the exponent a to vary between 1 and 2, one obtains a total mass in the Jupiter-Saturn (J-S) 
region in excess of 40 Mq, and as high as 180 M®. Therefore, it is likely that the region originally 
contained significantly more solids than ended up in the cores of Jupiter and Saturn. 

3. Implications of oligarchic growth 

Runaway growth (eg. Wetherill &; Stewart 1989, Kokubo & Ida 1996) transitions to a slower, 
self-limiting mode called "oligarchic growth" (Ida & Makino 1993, Kokubo & Ida 1998, 2000) when 
the largest protoplanets are still orders of magnitude less than an Earth mass everywhere in the 
nebula. Oligarchic growth has previously only been demonstrated to take place interior to about 3 
AU (Weidenschilling &; Davis 2000). Though Kokubo and Ida make estimates of protoplanet mass 
and growth timescales in the giant planet region, they point out that their simulations are restricted 
to annuli which are narrow compared to their radii, and thus cannot make strong predictions 
about how oligarchic growth works over a wide range in semimajor axis (Kokubo and Ida 2000). 
However, their most recent simulations (Kokubo & Ida 2000b) span a range of 0.5 to 1.5 AU, and 
show oligarchic growth proceeding in an outward-expanding "wave" over time. Also, Thommes 
(2000) and Thommes, Duncan &; Levison (2001) perform numerical simulations which suggest that 
oligarchic growth does take place in the outer Solar System, and proceeds on approximately the 
timescales predicted using the approach of Kokubo & Ida (2000). 
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Assuming that protoplanets remain approximately evenly spaced in Hill radii while growing — 
as is characteristic of oligarchic growth — their final mass is given by 



where a is the protoplanet semimajor axis, n is the spacing between adjacent protoplanets in Hill 
radii vh, and p is the fraction of the total mass in the zone [a — nrij/2, a + nri7/2] incorporated into 
the protoplanet (Kokubo and Ida 2000). Using numerical simulations, Kokubo & Ida (1998) show 
that n ~ 5 — 10. Models of giant planet formation by concurrent planetesimal and gas accretion 
suggest a surface density profile a oc a^^ (Pollack et et 1996). Using this in Eq. 3, one obtains 
a protoplanet mass independent of semimajor axis. Adopting uo = 10 g/ cm? and a spacing of 
n=7.5 Hill radii, one must then set the accreted mass fraction to about 0.8 in order to obtain a 
protoplanet mass of 10 . With this spacing, between three and five such bodies fit between 5 
and 10 AU. It is likely, therefore, that this region originally contained more than just the future 
solid cores of Jupiter and Saturn. At the same time, the more recent simulations of Thommes, 
Duncan and Levison (2001, preprint) indicate that accretion is quite inefficient in the outer Solar 
System; even in the Jupiter-Saturn region, the above value of p is probably overly optimistic, and 
a correspondingly higher density of solids in this region may have been necessary to form gas giant 
cores. But as we find (Section 5.5), a higher-density disk actually tends to increase the "success 
rate" of the model. 



In the nucleated instability picture of gas giant formation, there are three distinct phases of 
growth (Pollack ct al 1996): In Phase 1 (which itself encompasses the sub-phases of runaway and 
oligarchic growth; see Section 3 above), a solid core grows until it depletes most of the material 
in its feeding zone. Phase 2 is characterized by much slower growth, with gas accretion gradually 
coming to dominate over solids accretion. After a time of order several Myrs, the protoplanet 
contains comparable masses of gas and solids; around this time, the third phase of runaway gas 
accretion sets in and proceeds on a timescale of ~ 10^ years. One of the protoplanets in the Jupiter- 
Saturn region must have been the first to reach this point. Shorter formation timescales at smaller 
heliocentric radii argue for Jupiter having formed before Saturn. The relatively long time spent 
by a protoplanet on the "plateau" of Phase 2 means that several protoplanets could plausibly find 
themselves there by the time the first of them makes it to Phase 3. Furthermore, if the winner has 
a broad enough margin of victory, its rivals will have only had time to accumulate a few of gas, 
thus resembling the present-day ice giants in both solids and gas mass. 

A body increasing its mass from 10 M® to Jupiter's mass, Mj = 314 M®, expands its Hill 
radius by a factor of about three. The adjacent protoplanets will therefore have a high likelihood of 
being gravitationally scattered. In a three-body system — Sun, Jupiter and a single protoplanet — a 




(3) 



4. Scattering and circularization of the protoplanets 
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Jupiter-crossing protoplanet would continue to have close encounters with its much more massive 
scattcrcr, and so would remain coupled to it, unless it were scattered onto an unbound orbit and left 
the system altogether. In reality, however, one expects that as it crosses beyond the Jupiter-Saturn 
region, a body will encounter a less accretionally evolved part of the protoplanetary disk, consisting 
predominantly of much smaller bodies. As a result, the protoplanet will experience dynamical 
friction. This will tend to reduce the eccentricity of its orbit. If the eccentricity decays enough, 
the protoplanet's perihelion will be lifted away from Jupiter, thus decoupling it from its scattcrcr. 
Insofar as one can neglect interactions with other scattered protoplanets, its eccentricity will then 
monotonically decrease until it reaches equilibrium with the planctesimals. The eccentricity rate 
of change due to dynamical friction on a body of mass M in a swarm of mass m bodies can be 
expressed as (Weidenschilling et al 1997) 

^ = C{m{el) - Mel,)K, (4) 

and a similar expression gives the time evolution of the inclination, (e^)^/^ is the RMS eccentricity 
of the mass m bodies, and Kg is a definite integral depending on the inclinations and eccentricities, 
defined in Stewart & Wether ill (1988). The coefficient C is given by 

4((e?) + (ei))3/2' 

where G is the gravitational constant, vk the local Keplerian velocity, Pswi the spatial mass density 
of the swarm, and A is approximately the ratio of the maximum encounter distance between the 
body M and a member of the population of the swarm, to the maximum separation that results in 
a physical encounter (Stewart & Wether ill 1988). 

An equilibrium is therefore reached when 

ei,M ~ {el)m (6) 

This also means that as long as 

ei^M » (e^)m, (7) 

the eccentricity decay rate of the body M will only depend on the spatial mass density of the 
planetesimal disk, and will be essentially independent of the individual masses of the planetesimals. 
This feature will be exploited in the numerical simulations below. 



5. N-body simulations 

Putting the pieces together, one can now envision a scenario in which a) a rapidly growing 
Jupiter scatters its smaller neigbours outward, and b) these "failed cores" are decoupled from 
Jupiter and ultimately evolve onto circular, low-inclination orbits in the outer Solar System. TDL99 
performed three series of numerical simulations to test this model; about half of the simulations 



-6- 



produced a stable final configuration qualitatively similar to the present-day Solar System. Here, 
we expand on the earlier work, performing simulations to assess the effect of other planetesimal disk 
density profiles and of proto-Jupiter being initially not the innermost giant protoplanct. We also 
investigate how sensitive this scenario is to the timing of Saturn's final stage of growth, relative 
to that of Jupiter. Simulations are once again performed using the proven SyMBA symplectic 
integrator (Duncan, Levison & Lee 1998). This integrator is able to handle close encounters among 
massive bodies while preserving the symplectic properties of the method of Wisdom &; Holman 
(1991). 



5.1. Initial conditions: Set 1 

For the baseline set of simulations, a planetesimal disk of surface density 

(71 = 10(a/5 AU)-"^ g/cm^ (8) 

is used, with the disk extending from 5 to 60 AU. The above density profile is steeper than those 
of TDL99 (cr cx a'^ and a^^/^). The total disk mass between 10 and 60 AU (106 M®) is about half 
that of the first two sets of runs in our earlier work (216 M®), and similar to that in the third run 
(119 M0). The truncation at 60 AU is to keep the number of bodies in the simulations tractably 
small; in reality the planetesimal disk may have extended for hundreds of AU, far beyond the 
presently visible Kuiper belt. Observational evidence does point to the possibility of a truncation 
of the belt at ~ 50 AU (eg. Chiang &; Brown 1999). But if a "Kuiper cliff" exists, it is not 
necessarily primordial. 

Since we truncate the disk at a location where the planetesimal surface density is already very 
low (0.07 g/cm2 for the density profile of Eq. 8), one expects this not to have a strong effect on the 
evolution of any scattered protoplanets which cross the region r > 60 AU. As a test, we perform 
four pairs of runs comparing the evolution of a 10 body with initial o = 100 AU, e = 0.9 — 0.95 
in planetesimal disks with a surface density as prescribed by Eq. 8, in one case extending to 60 AU, 
in the other to 200 AU. No systematic difference in the 10 body's semimajor axis evolution 
over the first few Myrs is apparent between the two different disk sizes. The eccentricity, however, 
tends to be damped more rapidly in a 200 AU disk. We expect, therefore, that the subsequent 
runs somewhat underestimate the effectiveness with which those scattered protoplanets which cross 
beyond 60 AU are circularized, if the disk was in fact larger. In particular, a planetesimal disk with 
a radius of hundreds of AU may allow the retention of protoplanets which in our runs are scattered 
strongly enough to become unbound; we will explore this possibility in future work. 

The simulated disk is made up of equal-mass "planetesimals" , each having a mass of 0.2 M®. 
At twice the mass of Mars, these bodies far exceed the actual characteristic mass of planetesimals 
in the early outer Solar System. In reality planetesimals are thought to have been on the order of 1 
to 100 km in size, thus with a mass of ~ 10~^^ — 10~^ (eg. Lissauer 1987). The unrealistically 
large masses are chosen to keep the number of bodies manageably low, at slightly over 500. As 
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Fig. 1. - Eccentricities (top) and inclinations (bottom) in the outer Solar System at the present 
epoch, showing the giant planets as well as all Kuiper belt objects and Centaurs (objects with a 
< 30 AU) which have been observed at multiple oppositions, as of September 2001. The curve in 
the top panel shows the locus of orbits with perihelia at the semimajor axis of Neptune. KBO and 
Centaur data is taken from the Minor Planet Center site, cfa-www.harvard.edu/iau/mpc.html 
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mentioned in Section 4, when the large bodies have eccentricities high enough that Eq. 7 is satisfied, 
the eccentricity decay timescale will be effectively independent of the planetesimal masses. Thus 
despite the large planetesimal masses, the initial strength of eccentricity damping will be realistic. 
The equilibrium eccentricity condition Eq. 6, however, does depend on the planetesimal mass, so the 
equilibrium eccentricity reached by a large body among the planetesimals in the simulation will tend 
to be unrealistically high. Of course, a true equilibrium between protoplanets and planetesimals 
will not exist anyway, since mutual perturbations among the protoplanets will also have an effect. 

The initial planetesimal eccentricities and inclinations are given a Rayleigh distribution in e 
and i, (Kokubo and Ida 1992) with (e2)V2 = o.05, (i2)V2 = o.025 = 1.4°. In the numerical 
integrations, the planetesimals are treated as a "second-class", non-self-interacting population. 
Thus they are perturbed in their Keplerian orbits only by forces from the protoplanets, not each 
other. The protoplanets, on the other hand, are subject to forces from each other as well as from 
the planetesimals. This serves two purposes: It makes the simulations run much faster, since 
for N second-class bodies, the computation time scales as N instead of iV^. Also, it prevents 
unrealistically strong self-stirring of the disk. Of course, not modeling planetesimal interactions 
means that collective planetesimal effects are not accounted for. Wave phenomena could have had 
an important effect on the evolution of the planetesimal disk velocity distribution (eg. Ward and 
Hahn 1998), provided the disk was sufficiently massive and dynamically cold. However, it is unlikely 
that significant wave phenomena could persist once the initial scattering has taken place and the 
planetesimal disk has been stirred by eccentric 10 bodies. 

Also unmodeled is nebular gas, either as a source of aerodynamic drag (relevant for small 
planetesimals; eg. Adachi, Hayashi and Nakazawa 1976), or as a source of tidal forces (relevant for 
bodies > 0.1 M^; eg. Ward 1997). The former mechanism will keep the planetesimal disk more 
dynamically cold, but this does not affect the strength of dynamical friction on larger bodies until 
e^M ~ e^m. The latter effect is thought to cause the inward migration of 1 - 10 objects 
on timescales short compared to their formation times. This of course constitutes a potentially 
severe problem not only for our model, but for any model of giant planet formation which requires 
the accumulation of large solid cores. Addressing this problem is beyond the scope of our present 
work, but we summarize some reasons why it may in reality have been less severe, in Section 7. 

The 10 bodies are given a density p = 0.25/9® = 1.5 g/cm^, roughly equal to that of Uranus 
and Neptune. This gives them radii of 2.18 x 10^ km. Four such bodies are put in the simulation, 
initially on nearly circular and uninclined orbits. The orbits are spaced by 7.5 mutual Hill radii, 
starting from an innermost distance of 6 AU to allow for later inward migration by Jupiter. Thus 
the bodies' initial semimajor axes are 6.0 AU, 7.4 AU, 9.0 AU and 11.1 AU. Between 5 and 12 AU, 
the disk is depleted in planetesimals so that the surface density is still given by Eq. 8. Since the 
large bodies are spaced proportionally to their semimajor axes, their distribution is consistent with 
a surface density oc a~^. 

Integrating the surface density from 5 to 12 AU, the total mass is 51.5 M®. With 40 M® of this 
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Fig. 2. — Initial state for runs in Set 1, showing eccentricity (top) and inclination (bottom) versus 
semimajor axis. The larger circles denote the four 10 protoplanets, and each of the small dots 
represents a 0.2 planetesimal. The planetesimal density in the vicinity of the protoplanets is 
decreased to keep the density of protoplanets plus planetesimals consistent with the surface density 
given by Eq. 8. 
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in the large bodies, this leaves 10.5 in planetesimals in the Jupiter-Saturn region, in addition to 
94.7 between 12 AU, and the outer disk edge at 60 AU. In summary, then, the initial conditions 
amount to a state where Phase 2 of giant planet formation has been reached between 5 and 12 
AU, with ~ 80% of the planetesimals having accreted through oligarchic growth into four bodies 
of 10 each, while no large bodies have yet formed beyond this region (Fig. 2). The adoption of 
equal-mass protoplanets is a simplification; apart from stochastic variations among the oligarchic 
growth endproducts, one can expect some intermediate bodies, with masses perhaps up to a few 
M®. Such bodies are likely to ultimately be cleared, along with the planetesimals, from the giant 
planet region. However, they may end up playing a role in the dynamics of the trans-Neptunian 
region; see Petit, Morbidelli and Valsecchi (1999). 

For computational reasons, the inner simulation radius is chosen as 1 AU; any body whose 
orbit penetrates this boundary is eliminated from the system. This is done because a limitation 
of the SyMBA integrator used here is its inability to handle close perihelion passages. Although 
a new version of SyMBA has since been developed which removes this restriction (Levison and 
Duncan 2000), the runs presented here predate this development. In any case, this limitation has 
little relevance for the runs presented here; typically, less than ten planetesimals are lost at the 
inner boundary over the course of a 5 Myr simulation. The base timestep is chosen as 0.05 years, 
giving 20 steps per orbital period for an orbit with its semimajor axis at the inner radius, and over 
200 steps per orbit for Jupiter. Experimentation shows that this timestep is small enough that the 
energy of the system is well-conserved. Runs initially go to 5 Myrs; in cases where the system still 
appears to be undergoing rapid evolution at this point, the runs are extended by another 5 Myrs. 



5.2. Set 1 results 

To model gas accretion, SyMBA was modified to allow a subset of bodies to have artificially 
time-varying masses (in addition to any changes in mass resulting from the accretion of other 
bodies). For the runs in Subset 1, it is assumed that the innermost protoplanet undergoes runaway 
gas accretion first, and grows into Jupiter. The simulations start at the time when this happens, 
which should be a few million years into the life of the solar system (see above). Runaway gas 
accretion is simulated by increasing the body's mass over the first 10^ years of simulation time, 
from its original mass of 10 M0 to 314 M0, approximately the present Jupiter mass. A linear 
growth in mass is used; this is deemed appropriate since the actual time evolution of mass during 
runaway growth is highly uncertain. Also, as the simulations will show, 10^ years is roughly the 
response time of the system, and the system's subsequent evolution is therefore unlikely to be 
affected in a systematic way by the exact form of the time evolution of the runaway-phase mass 
growth. 

Set 1 consists of eight alternate realizations of a run, differing only in the initial phases of the 
four protoplanets; for each, the angles J7, uj and M are randomly generated. As will be seen, the 
stochasticity of the system ensures that this difference in phases is sufficient to bring about a very 
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different evolution in each of the versions of the run. 

Fig. 3 shows the evolution of one of the eight runs, denoted as Run IF, which after 5 Myrs 
produces final protoplanet orbits that bear a particularly close resemblance to those of the giant 
planets in the present-day Solar System. Semimajor axis versus time is plotted for each of the four 
protoplanets; the innermost one has grown into Jupiter after the first 10^ years. By this time, the 
protoplanet orbits begin to mutually cross, and strong scattering occurs. The protoplanet plotted 
in red briefly has its semimajor axis increased to greater than 100 AU. However, dynamical friction 
acts to reduce eccentricities universally, decoupling the protoplanets from Jupiter and from each 
other. By about 1.2 Myrs, none of the protoplanets are on crossing orbits anymore. After about 
3 Myrs, the bodies no longer undergo any changes in semimajor axis greater than a few AU on a 
million-year timescale. At this point the orbits are well spaced and all eccentricities are < 0.05, 
with no large fluctuations. Due to the large stochastic variations among runs, the exact timescales 
differ, but we find that orbits typically become noncrossing after less than 5 Myrs in these and 
subsequent runs. 

Subsequent semimajor axis evolution proceeds by scattering of planetesimals by protoplanets, 
rather than scattering of the protoplanets off each other and Jupiter. As planetesimals are scattered 
among Jupiter and the protoplanets, the former experiences a net loss of angular momentum while 
the latter experience a gain. Thus Jupiter's orbit shrinks, while those of the protoplanets expand 
(Fernandez and Ip 1996, Hahn and Malhotra 1999). This phase takes place over a timescale of 
several tens of Myrs and ends, at the very latest, when the planetesimals have been cleared from 
among the planets. Because the length of migration during this phase is only a few AU, and to save 
time, most of the runs are stopped after 5 Myrs. As an example. Fig. 4 shows Run IF continued 
to 50 Myrs; the net migration of the outer two protoplanets subsequent to 5 Myrs is only ~ 1 - 4 
AU outward. 

The semimajor axes at which the scattered bodies end up are very noteworthy, if one compares 
them to the present orbits of the giant planets (Fig. 1). At 5 Myrs, the outer two protoplanets are 
at 16 and 26 AU. The innermost one is at 11 AU, while Jupiter is at 5.5 AU. This configuration of 
orbits is very similar to that of the present Solar System, where Uranus and Neptune are at 19 and 
30 AU respectively, Saturn is at 9.6 AU, and Jupiter is at 5.2 AU. And at 5 x 10^ years, after some 
more net outward migration, the outer two protoplanets' semimajor axes are even closer to those 
of Uranus and Neptune (though "proto-Saturn" , having also moved outward, is further away from 
its present orbit). Eccentricities and inclinations are likewise very close to their present values. 

The end states of all the runs in Set 1 are summarized in Fig. 5. Depicted are snapshots of 
eccentricity versus semimajor axis at 5 Myrs, except for Runs lA and IH, which were continued 
to 10 Myrs. These were extended because one or more of the protoplanets still had a high but 

decreasing eccentricity at 5 Myrs. 

Runs lA, ID, IF, IG and IH result in a final ordering of orbits that is at the very least broadly 
consistent with the present Solar System: Jupiter is the innermost body, with the other three bodies 
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Fig. 3. — Run IF: Evolution of semimajor axis (bold lines), perihelion distance q (thin lines) and 
aphelion distance Q (dotted lines) of the four 10 Mq protoplanets. The protoplanet which grows 
to Jupiter mass (314 Mq) over the first 10^ years of simulation time is shown in black. 



- 13 - 




Fig. 4. — Run IF continued to 50 Myrs. Between 5 Myrs and 50 Myrs, the net migration for 
Jupiter and the protoplanets, going from inside to outside in semimajor axis, is -0.2 AU, 1.5 AU, 
4 AU, and 1.3 AU. 
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Fig. 5. — Endstates of the eight Set 1 runs, after 5 Myrs of simulation time, except for lA and 
IH, which were continued on to 10 Myrs. Eccentricity is plotted versus semimajor axis. The three 
different sizes of points denote planetesimals (smallest), 10 protoplanets (medium), and Jupiter 
(largest). Planetesimal orbits crossing Jupiter or any of the protoplanets arc generally unstable on 
timescales short compared to the age of the Solar System, thus the region among the protoplanets 
would be essentially cleared of planetesimals long before the present epoch. 
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interior to the region of the Kuiper belt, and eccentricities low enough that no protoplanets cross 
each other or Jupiter. Out of these five runs, ID and IF in particular resemble the present Solar 
System. Of course, to actually reproduce the Solar System, another important event has to take 
place: The next protoplanet beyond Jupiter must also undergo a runaway gas accretion phase to 
acquire an envelope of mass > 80 and become Saturn. With the innermost core on a stable 
orbit in the vicinity of Saturn's present location, though — as it is in runs ID and IF — the time 
delay between Jupiter and Saturn's runaway phases is not strongly constrained. A more detailed 
investigation of the role of Saturn's growth will follow in Section 5.7. 

In contrast, Uranus and Neptune must somehow have been prevented from later reaching 
runaway gas accretion. One possibility is that they simply ran out of time, still caught in the 
long plateau of the second giant planet growth phase (Pollack et al 1996; see Section 4) when the 
nebular gas was removed. Alternatively, it may be that the gas disk was truncated early on by 
photoevaporation between the orbit of Saturn, and the eventual orbit of Uranus (Shu, Johnstone 
and HoUenbach 1993). This possibility is discussed further in Section 7. 

The final state of the planetesimal disk differs substantially among the runs. In IB, IC and 
IE, the planetesimal disk is largely unperturbed over most of its radial extent, while in IF and 
IH, eccentricities have been greatly increased throughout the entire disk. lA, ID and IG are 
intermediate cases. The extent of the perturbation simply depends on how much of the disk is 
crossed by the protoplanets over the course of the run. In Run IF, for example. Fig. 3 shows that 
one body's aphelion spends some time beyond 100 AU. On the other hand, in IC, no protoplanet's 
aphelion ever goes further out than Qmax — 18 AU. In those runs where Qmax < 60 AU (the disk 
radius), Qmax corresponds closely to where the planetesimal disk makes a transition from perturbed 
to unperturbed. The outer limit of the planetesimal disk's excitation by scattered protoplanets has 
been referred to as the "fossilized scattered disk" (TDL99). The contemporary scattered disk, 
by contrast, consists of objects which were scattered by Neptune after the latter had attained its 
current orbit (Duncan and Levison 1997). A more detailed discussion of the fossilized scattered 
disk follows in Section 6.2. 

Three of the runs produce systems at 5 Myrs that are irreconcilably different from our own. In 
IB, one of the protoplanets has merged with Jupiter. In IC, a protoplanet has been scattered onto 
an orbit interior to Jupiter, in the region of the present-day asteroid belt, with its semimajor axis 
at 3.4 AU, its perihelion at 2.8 AU and its aphelion at 3.9 AU. It attains a stable orbit not crossing 
Jupiter even though the only planetesimals available for damping interior to Jupiter arc those few 
which arc also scattered there. This is possible because the protoplanet is scattered not just by 
Jupiter, but also by the other protoplanets. Even assuming the protoplanet could be subsequently 
removed from this region, such an event would most likely have cleared much of the asteroid belt. 
Finally, in IE, two of the protoplanets have merged. It should be noted that mergers — or, indeed, 
ejections — which reduce the number of protoplanets are not intrinsically a problem, since extra 
ones may have existed. Scenarios with five protoplanets (Jupiter + Saturn + 3) will be explored 
in Section 5.4. 
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5.3. Set 2: Dependence on initial ordering 

How strongly does the final configTiration of the system depend on the initial ordering of the 
protoplanets? The next set of runs, Set 2, uses the same initial conditions, within a random 
variation in the protoplanet phases, as Set 1. However, for these runs it is the second-innermost 
protoplanet, rather than the innermost one, which undergoes simulated runaway gas accretion. 
One can reasonably expect that this will favour an outcome like IC (Fig. 5), where a protoplanet 
is scattered inward into the region of the asteroid belt. 

The end states of the runs are shown in Fig. 6. In six of the eight runs, a protoplanet has 
indeed ended up interior to Jupiter. However, in two cases (2C and 2D), Jupiter is the innermost 
body. Thus it appears that if Jupiter does not grow from the innermost protoplanet, the likelihood 
of ending up with a final configuration similar to our Solar System declines, though such an outcome 
continues to be quite possible. 

5.4. Set 3: Dependence on number of cores 

How sensitively does the end state depend on the initial number of core-sized bodies? In 
the next set of simulations, an extra 10 protoplanet is added. All protoplanets are more 
tightly spaced, by 6.5 instead of 7.5 mutual Hill radii. Starting, again, from 6.0 AU, the outermost 
protoplanet is therefore initially at 12.2 AU. The surface density of planetesimals in the region of 
the protoplanets is reduced to keep the average surface density unchanged at ai = 10(a/5AU)~^ 
g/cm^. The innermost protoplanet is, again, increased in mass to that of Jupiter over the first 10^ 
years. 

Fig. 7 shows the endstates of the runs. This time all runs initially have a length of 10^ years, 
since the larger number of bodies take longer to decouple from each other. Run 3A is continued 
to 1.5 X 10^ years, because after 10^ years some of the protoplanets are still on crossing orbits. 

Eccentricities are uniformly low and the protoplanet orbits arc well-spaced for the most part, thus 
the systems have a good chance of being stable indefinitely. However, once all the planetesimals 
have been scattered from among the protoplanets, so that the latter are no longer subject to 
dissipative forces, some of these systems may still become unstable. This caveat applies to all the 
runs presented in this work, but generally speaking, larger numbers of bodies increase the potential 
for instability (Levison, Lissauer and Duncan 1998). 

All of the protoplanets remain in six of the eight runs, thus resulting in systems with one 
too many planets relative to the present Solar System. However, in Runs 3B and 3F, one of the 
protoplanets is ejected from the Solar System, leaving the right number of bodies behind. In both 
of these runs, a protoplanet ends up with a semimajor axis within 10% of present-day Uranus and 
Neptune, respectively, though both "Saturns" are too far out. One may conclude that with one 
extra initial protoplanet in the Jupiter-Saturn region, scattered protoplanets continue to be readily 
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circularized, and the resulting systems tend to look like ours with one extra outer planet. However, 

a system with four giant planets remains a possible outcome. Also, the subsequent formation of 
Saturn from one of the protoplancts may trigger more ejections, especially in those cases where the 
inner protoplanets are more closely spaced, such as 3D and 3H. 



5.5. Set 4: A more massive planetesimal disk 

In this set of runs, a number of parameters are changed to simulate a system with a more 
massive planetesimals disk. The protoplanets are now 15 bodies. The planetesimal disk surface 
density profile is still cx a~^, but it is scaled up to be 15 g/cm^ at 5 AU; that is, 

0-2 = 15(a/5 AU)~2 (9) 

Other minor differences are a legacy of chronologically earlier runs. The innermost protoplanet is 
initially at 5.3 AU, and successive protoplanets are spaced by only 5.8 mutual Hill radii. Thus the 
outermost protoplanet is initially at 9.0 AU. The individual planetesimals have a mass of 0.24 Mq. 

The innermost protoplanet, as before, has its mass increased to 314 over the first 10^ 
simulation years. The end states of the runs are shown in Fig. 8. All except C yield the correct 
number and ordering of bodies, and eccentricities are uniformly low. This "success rate" is higher 
than that of Set 1, in which only five out of eight runs yield qualitatively the correct orbital 
configuration. This is accounted for by the more massive planetesimal disk; it provides stronger 
dynamical friction, so that scatterings of protoplanets tend to be less violent, and subsequently, 
orbits tend to be circularized and mutually decoupled more quickly. Runs A, D and H end up with 
protoplanet orbits that are particularly close to those of Saturn, Uranus and Neptune. Jupiter 
systematically ends up at too small a heliocentric distance, indicating that the initial distance of 
5.3 AU is too small. Also, the larger disk mass gives the protoplanets and Jupiter more planetesimals 
to scatter, and thus increases the distance they travel due to angular momentum exchange (Section 
5.1). 



5.6. Set 5: A shallower disk density profile 

In this set of runs, a shallower planetesimal disk surface density, oc a~^/^, is used. The disk 
now begins at 10 AU, with the interior region being initially occupied solely by the protoplanets. 
In other words, it is assumed that in the epoch from which the runs start, all but a negligible mass 
of the planetesimals among the protoplanets has been swept up or scattered from the region. The 
surface density profile is given by 

(73 = 1.8(a/10 AU)-^/^ g/cm^ (10) 

and thus a total mass of 123 is contained in the disk between 10 and 60 AU. Extending this 
planetesimal disk inward to 5 AU would yield a surface density there of only 5 g/cm^, and a total 
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mass between 5 and 10 AU of only 25 M®. Despite this, we still put four 15 in this region; 

thus, it is assumed that the original planetesimal surface density profile in this region was steeper, 
perhaps due in part to redistribution of water vapor from other parts of the disk to the vicinity of 
the snow line (Stevenson & Lunine 1988). 

This set consists of twelve runs, each to 5 Myrs. As before, the inner body's mass is increased 
to 314 Mm over the first 10^ years of simulation time. The endstates are shown in Fig. 9. Eight out 
of twelve runs possess the right number and ordering of bodies. Eccentricities of the protoplanets 
and Jupiter are < 0.1 in five of these. This "success rate" is close to that of Set 1 (which has 
a similar total disk mass); this model thus does not appear to be highly sensitive to changes in 
density profile alone. 

Various degrees of disruption of the planetesimal disk can again be seen. Most show a sharp 
transition between a disrupted region crossed by the scattered planetesimals, and a largely undis- 
turbed outer region. In Run 5A, for example, this transition occurs at slightly below 40AU, while 
in Run 5F, it is located between 45 and 50 AU. Runs 5B and 5L show strong disruption throughout 
the entire disk, indicating that all of it was crossed by one or more protoplanets. Another state can 
be seen in Run 5G, where all eccentricities in the outer part of the belt are uniformly raised. This 
occurs when a protoplanet crosses the outer disk with an inclination high enough that it spends 
most of its orbit above or below, rather than inside, the disk. The disk planetesimals are then 
excited primarily by long-range secular effects rather than by short-range scattering encounters 
(TDL99). Another example of this effect can be seen in Run IG above. 

5.7. Set 6: The role of Saturn 

Thus far, the only gas giant in the simulations has been Jupiter. The innermost of the scattered 
protoplanets does tend to end up near the present location of Saturn. However, to reproduce the 
Solar System, it must at some point accrete ~ 80 Mq of nebular gas. We therefore investigate what 
effect the subsequent growth of a Saturn-mass object has on our model. 

A Set 1 run which produced a good Solar System analogue (IF; Fig. 3) is used as a starting 
point. At about 1 Myr, the protoplanets and Jupiter are no longer on crossing orbits. At this 
time, the initially outermost protoplanet (plotted in blue) has become the innermost one, closest 
to Jupiter at 10 AU. We perform a set of five runs which branch off from this point. In these 
runs, the innermost protoplanet has its mass increased to that of Saturn over a 10^ year interval, 
starting at 1, 1.2, 1.4, 1.5 and 1.6 Myrs, respectively. We choose to start only after the protoplanets 
are on noncrossing orbits in order to avoid ending up with an eccentric "Saturn" ; the eccentricity 
evolution of an initially eccentric giant planet in a gas disk is uncertain (Lin et al 2000). 

The endstates are shown in Fig. 10. No protoplanets have been lost from the system by the 
end of the runs. A protoplanet still has a high eccentricity in the 1 Myr case; this is because 
this protoplanet's perihelion is still very close to the innermost protoplanet's orbit at 1 Myr, and 
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it suffers strong perturbations as the latter grows to Saturn's mass. In the other cases, however, 
Saturn's formation does not cause large eccentricities in the protoplanets. This is as one would 
expect; at its final mass, and at 10 AU, Saturn's Hill radius is 0.45 AU, and by 1.2 x lO*' years, 
the closest protoplanet's perihelion is at ~ 14 AU, almost 9rH away, thus unlikely to be in reach 
of strong scattering. We leave for future work the effect of Saturn's gas accretion on systems with 
more protoplanets, such as those in Set 3 in Section 5.4 above. Such systems tend to have weaker 
stability, and thus may be more susceptible to disruption by Saturn's final growth spurt. 

One effect visible in Fig. 10 is that the semimajor axis of Saturn at 5 Myrs tends to be smaller 
than that of the innermost scattered protoplanet — the putative proto-Saturn — in those runs where 
the gas accretion of Saturn is not modeled (see for example Fig. 5). When a protoplanet grows to 
Saturn's mass, its subsequent migration speed is much slower, since the rate of migration depends 
on how much mass in planetesimals it scatters relative to its own mass. This counteracts the 
tendency of the innermost scattered protoplanet to end up at a semimajor axis larger than that of 
present-day Saturn in the other runs, where only Jupiter grows. 

6. Scattered protoplanets and the small body belts 

The simulations presented above show that scattering of giant planet core-sized protoplanets 
is a violent event, which leaves a strong dynamical signature on the surrounding planetesimal disk. 
The asteroid belt and the Kuipcr belt are therefore the natural places to look for evidence of large 
scattering events in the Solar System's early history. 

6.1. The asteroid belt 

It is those members of the asteroid belt larger than about 50 km in diameter which are of 
interest in inferring properties of the early Solar System; smaller bodies cannot be primordial 
because they could not have survived intact for the age of the Solar System (eg. Petit, Morbidelli 
and Valsecchi 1999). This population displays a degree of dynamical excitation that is not readily 
explained. The other puzzling feature of the asteroid belt is its severe mass depletion — by at the very 
least a factor of 10^ — relative to the amount of mass the region originally contained, extrapolated 
from the terrestrial zone (Weidenschilling 1977). 

In the context of the model presented here for the early evolution of Uranus and Neptune, 
the initial violent scattering of protoplanets is perhaps the most obvious candidate to look to for 
perturbation of the asteroid belt. Indeed, in numerous runs, protoplanets spend some time interior 
to the orbit of Jupiter, crossing part of the asteroid belt region for up to ~ 10^ years at a time. 
One might expect that such an occurence would wreak havoc in the the asteroid belt. To explore 
this possibility, we perform six additional runs with planetesimals added in the asteroid belt region. 
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Fig. 11. — Eccentricities and inclinations of planetesimals in the asteroid belt region, interior to 
Jupiter (the large dot), at 3 x 10*^ years (top panel) and 2 x 10""' years (bottom panel). These 
times are, respectively, just before and just after the period during which a protoplanet repeatedly 
crossed the region interior to Jupiter. Jupiter in this run has moved inward to ~ 4.8 AU, 0.4 AU 
less than its present semimajor axis. The curve marks the locus of Jupiter-crossing orbits 
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The individual bodies have a mass of 0.024 Mq, and are distributed with a surface density of 

abeit = 8.0(a/l AU)-i g/cm^ (11) 

between 2.5 and 4.5 AU. This shallow density profile is the same as the one used by Chambers and 
Wetherill (1998) in the terrestrial region, which in turn was chosen to be more consistent with the 
large densities required at larger heliocentric distances to form Jupiter and Saturn. The protoplanet 
masses in this case are 15 M®. We find that eccentricities can get excited to their present values 
in this way, though only down to the crossing protoplanet's minimum perihelion distance, which 
in none of the runs performed reaches the inner edge of the belt. Inclinations fare more poorly; 
protoplanets seldom raise them much above 10°, which is less than the present median inclination 
of asteroids beyond 2.5 AU. Also, very little mass is scattered out of the belt while the protoplanets 
are crossing it. Fig. 11 shows "before and after" snapshots of the run which, out of the six, displays 
the strongest disruption of the region interior to Jupiter by scattered protoplanets. As can be seen, 
few bodies, apart from those that are nearly Jupiter-crossing, attain inclinations above 15°. 

Scattered protoplanets may have had a more indirect role in the excitation and mass depletion 
of the asteroid belt. If the majority of the solids in the asteroid belt accreted to form a system 
of planet-sized bodies (raising eccentricities and inclinations of the remaining planetesimals in the 
process), scattered protoplanets crossing into the region may have contributed to making this 
system unstable, akin to the model of Chambers (1999). Also, if proto-Saturn was still sufficiently 
close to Jupiter when it accreted its massive gas envelope, the subsequent migration of the gas 
giants may have been enough to sweep the inclination-exciting viq secular resonance through most 
of the asteroid belt (cf. Gomes 1997, Levison et al 2001) As an example, in one of the runs from 
Set 6, Saturn commences growing at 4 x 10^ years and reaches its final mass at 5 x 10^ years, at 
which point Jupiter and Saturn are at 5.7 and 8.4 AU, respectively. This places the i^iq resonance 
at ~ 3.5 AU, and it moves inward as Jupiter and Saturn move apart. In contrast, Gomes (1997), 
using a more moderate range of migration for Jupiter and Saturn, finds that only the region inward 
of 2.7 AU is crossed by the uie resonance. 



6.2. The Kuiper belt 

In the present Solar System, a new class of Kuiper belt object (KBO) has recently been 
identified (Duncan and Levison 1997, Luu et al 1997). These objects have semimajor axes and 
eccentricities such that they lie near the locus of Neptune-encountering objects shown in Fig. 1. 
They arc thought to be part of a population referred to collectively as the scattered disk — formerly 
low-eccentricity KBOs which have had their orbits changed by close encounters with Neptune. 
Many of the simulations in Section 5 show an analogous class of planetesimals in their "Kuiper 
belt" regions. However, these fall on the locus of orbits crossing not the final semimajor axis of the 
outermost protoplanet, but the furthest aphelion distance of any of the protoplanets during their 
initial high-eccentricity phase. Since these orbits are no longer being crossed by a protoplanet. 



they will be stable over long times. One can refer to these structures as "fossilized" scattered 
disks (TDL99), because they preserve part of the dynamical history of the planetesimal disk. Such 
structures only appear in runs where the initial scattering was strong enough that one or more 
protoplanets had their aphelia increased to well beyond the final semimajor axis of the (ultimately) 
outermost protoplanet. 

Observations of our Solar System's Kuiper belt do indeed reveal an anomalously high degree 
of excitation (eg. Petit, Morbidelli &: Valsecchi 1999, Malhotra, Duncan k, Levison 2000). The 
eccentricities and, to a lesser degree, inclinations of bodies in mean-motion resonances with Nep- 
tune, particularly the 2:3 resonance at 39.5 AU, can be explained by resonance sweeping during 
Neptune's migration, as can the paucity of objects on nonresonant orbits interior to 39 AU (Mal- 
hotra 1995). However, the high inclinations found beyond ~ 41 AU in what is commonly called 
the "classical" Kuiper belt, cannot be explained in this way. Petit, Morbidelli and Valsecchi (1999) 
propose large (up to 1 M®) Neptune-scattered planetesimals as the mechanism which stirred and 
cleared the belt. But even when such bodies remain in the belt for 100 Myrs, the inclinations they 
raise are almost always less than 20°. 

Can the excitation of the Kuiper belt be better accounted for if we are thus far only seeing 
the part of it that is interior to the locus of a fossilized scattered disk? The inclinations raised by 
a protoplanet can be directly obtained from the simulations of Section 5. The runs of Set 1, Set 3 
and Set 5 will be used for comparison. Inclinations for Sets 1, 3 and 5 are shown in Figs. 12, 13 
and 14, respectively. Inclinations beyond the outermost protoplanet are excited up to a maximum 
of about 30°, similar to those observed in the classical Kuiper belt today (see Fig. 1) However, 
observations of the Kuiper belt are biased against high-inclination objects. Brown (2001) derives 
a de-biased inclination distribution function for the classical belt: 



with a = 0.93 ± 0.02, cti = 2.2zb'|, and (T2 = 17 ± 3. One can define a parameter i' = cos~^(cos(i)) 
to give a measure of the characteristic inclination of a population. For the de-biased distribution 
above, i' = 21°. However, in the runs presented here, the largest i' in the region corresponding 
to the classical belt (between the outermost planet's 2:3 and 1:2 mean-motion resonance) is only 
15°. Thus, although higher inclinations are produced here than in the large Neptune-scattered 
planetesimals model of Petit, Morbidelli and Valsecchi (1999), the inferred full velocity distribution 
of the classical Kuiper belt still cannot be accounted for. 

It is possible to estimate with a simple numerical experiment if a planet as large as Uranus 
or Neptune can in principal excite the Kuiper belt to observed values. This experiment consists 
of a single Uranus-mass planet on an orbit with a = 45AU, e = 0.2, and i = 25°. The planet is 
embedded in a swarm of 400 masslcss test particles informally spread from 35AU to 55AU, with 
initial e = 0.01 and i = 1°. Fig. 15 plots the i' of the particles as a result of scattering off of 
the planet. It shows that a planet the mass of the ice giants can indeed excite the Kuiper belt to 
i' = 20° in a million years. 
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Fig. 12. — Counterpart to Fig. 5, showing inclination versus semimajor axis for the endstates of 
the runs in Set 1. 
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Fig. 13. — Counterpart to Fig. 7, showing inclination versus semimajor axis for the endstates of 
the runs in Set 3. 
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Fig. 14. — Counterpart to Fig. 9, showing inclination versus semimajor axis for the endstates of 
the runs in Set 5. 
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Fig. 15. — Characteristic inclination i' versus time of test particles between 35 and 55 AU, gravi- 
tationally stirred by a Uranus-mass body initially having a = 45 AU, e = 0.25, and i = 25°. The 
test particle inclinations reach the debiased value for the Kuiper belt, i' = 20° (Brown 2001), after 
about 1 Myr. 



-32- 



However, in none of the runs we performed was such a high inchnation imparted on a pro- 
toplanet. Alternatively, a less inclined protoplanet may be able to reproduce the Kuiper belt 
inclinations if it remains on a belt-crossing orbit significantly longer than 1 Myr (and thus longer 
than in our runs). This may require a less massive trans-Saturnian planetesimal disk, in order 
to increase the dynamical friction timescale, though that in turn would increase the chances of 
protoplanets being scattered out of the Solar System altogether. A longer circular ization time will 
also result if a protoplanet can be decoupled from Jupiter and the other protoplanets at a larger 
semimajor axis, so that the planetesimal surface density averaged over its orbit is lower. For the 
latter situation to have a better chance of occurring, the planetesimal disk needs to be extended to 
larger heliocentric distance. We will investigate this issue further in future work. 

7. Discussion 

The conventional picture of Uranus and Neptune's formation, whereby the ice giants accrete 
near their present heliocentric distances, has grave problems. Numerical simulations have not been 
able to produce ~ 10 M0 objects in the trans-Saturnian region in the lifetime of the Solar System 
without significantly increasing protoplanet radii, or invoking dissipational forces and planetesimal 
disk densities too large to be consistent with a physically plausible protostellar disk. 

Building on our previous work (TDL99), we have performed additional N-body simulations of 
the evolution of the outer Solar System starting at the time when the first gas giant forms. At 
this point, we assume that a number of 10 M0 objects have formed at a heliocentric distance of 
roughly 5 to 10 AU, as is suggested by the oligarchic growth model (Kokubo & Ida 1998, 2000). 
Using a variety of different initial conditions, we find as before that the accretion of Jupiter's gas 
envelope causes the remaining protoplanets to become violently unstable. In most cases they are 
scattered onto high-eccentricity orbits in the trans-Saturnian region. With most of its orbit now 
crossing the largely pristine trans-Saturnian planetesimal disk, a scattered protoplanet experiences 
dynamical friction and has its eccentricity rapidly damped. As a result, the protoplanets tend to 
end up on nearly circular, well-spaced orbits on a Myr timescale, with semimajor axes comparable 
to those of Saturn, Uranus and Neptune. Of the simulations which initially contain a total of 
four giant protoplanets and form Jupiter from the innermost, the majority produce final orbital 
configurations similar to that of our outer Solar System. Such systems are produced — though with 
lower probability — even if ones adds an extra protoplanet, or lets the second-innermost protoplanet 
grow into Jupiter. These results strengthen our earlier conclusion that if Uranus and Neptune 
shared the same birthplace as the gas giants, they could then readily have been delivered to their 
present orbits. 

The role of migration in the formation of Uranus and Neptune was previously investigated 
numerically by Ipatov (1991), based on an idea by Zharkov and Kozenko (1990). Ipatov concludes 
that planetary embryos of a few may have originated just outside the orbit of Saturn, to 
migrate outward and later grow into Uranus and Neptune, provided that they did not acquire 
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high eccentricities during this process. In contrast, we find that Uranus and Neptune could have 

originated from anywhere in the Jupiter-Saturn region, and that initially high eccentricities — which 
nearby bodies naturally tend to acquire during Jupiter's final growth phase — are in fact a powerful 
mechanism for rapidly transporting them outward. Also, the long growth timescales in the outer 
Solar System suggest that Uranus and Neptune likely already completed most of their growth in the 
Jupiter-Saturn region; even with a "head start" of a few Earth masses, the formation of Uranus- 
and Neptunc-niass objects much beyond 10 AU within the age of the Solar System is unlikely 
(Levison and Stewart 2001, Thommes, Duncan and Levison 2001) 

Can one find any observational support for this model in the present-day Solar System? The 
high inclinations in the classical Kuiper belt point to strong dynamical excitation in the past, 
and the simulations performed here do produce high inclinations in this region as a natural side 
effect. However, the simulations all fall short of reproducing the mean debiased inclination of 
the classical Kuiper belt. Strong observational support would be provided by the discovery of a 
fossilized scattered disk in the Kuiper belt, and a dynamically colder population beyond. It is 
tempting to link the trans-Neptunian object 2000 CR105 with a fossilized scattered disk; its high 
eccentricity (0.8) is characteristic of a scattered disk object, but recent observations (Gladman ct al 
2001) have established that its perihelion is at 44 AU, far beyond the reach of Neptune. However, 
none of the fossilized disk objects in our simulations acquire semimajor axes as high as that of 
2000 CR105 (216 AU). Fig. 16 shows a plot of perihelion distance versus semimajor axis for Set 
5, revealing only one case (B) where one or more objects simultaneously acquire a semimajor axis 
of ~ 100 AU and a perihelion distance significantly further out than the (circularized) outermost 
protoplanet. All other sets of runs fare even more poorly. 5B is a run in which a protoplanet spends 
a long time at high eccentricity, and excites particularly high planetcsimal inclinations in the disk. 
A long circularization time may thus be an important ingredient in reproducing both objects like 
2000 CR105 and the high inclinations of the Kuiper belt; we will address this possibility in future 
work. 

Findings regarding the deuterium to hydrogen (D/H) ratios of Uranus, Neptune and comets 
are particularly interesting in the context of this model. Prom infrared observations, Feuchtgruber 
et al (1999) find that the D/H ratios of the ice giants are lower than those of the Oort cIoTid comets 
Halley (Eberhardt et al 1995), Hyakutakc (Bockelee-Morvan et al. 1998) and Hale-Bopp (Meier et 
al 1998) by a factor of approximately three, a difference large compared to the uncertainties of the 
measurements. Oort cloud comets are thought to have originated primarily in what is today the 
Uranus-Neptune region (Duncan, Quinn k, Tremaine 1987, Fernandez 1997), and the notion of a 
common birthplace is supported by the comets' similar D/H ratios. Though a sample size of three 
is clearly very small, this discrepancy between the comets and the ice giants would seem to present a 
further problem for any scenario in which Uranus and Neptune form in the trans-Saturnian region, 
since they should then share the chemical composition of the comets. However, this is exactly what 
one would expect if the ice giants originally formed at a smaller heliocentric distance, where higher 
temperatures would have made for a lower D/H ratio. 
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Fig. 16. — Endstates of the Set 5 runs, showing perihehon distance versus semimajor axis. 
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An aspect not modeled in any of the simulations is the gravitational interaction of the bodies 
with a gaseous disk. It has been shown that gas disk tidal forces can cause rapid inward migration 
of protoplancts (eg. Ward 1997). In fact, the speed of migration may be peaked for objects of ~ 

10 M^, taking place on timcscales of 10^ years or less. This peak corresponds to the transition 
between so-called Type I migration, where a body's resonant interaction with the gas disk gives 
rise to a torque imbalance, to Type II migration, where the object opens a gap in the gas disk and 
is subequently locked to the disk's viscous evolution. Tidal migration therefore poses a problem 
for any model of giant planet formation: how do they avoid spiralling into the central star as they 
form? However, tidal torques may not in fact operate throughout the whole disk. Gammie (1996) 
develops an accretion disk model which, beyond ~ 0.1 AU, only transports angular momentum in 
a relatively thin surface layer. Thus the bulk of the disk would have zero viscosity, and objects 
embedded in it would be subject to neither Type I nor Type II orbital decay. Also, the disk may 
be truncated early on by photoevaporation from the central star (Shu, Johnstone and Hollenbach 
1993), or from surrounding stars, as is seen in the Orion Nebula proplyds (Johnstone, Hollenbach 
and Bally 1998). Bodies scattered beyond the truncation distance would then be safe from nebula 
tides, and those near the edge would experience a net positive torque and migrate outward. Of 
course, such a scenario has the added advantage that, as discussed in the above works, it accounts 
for Uranus and Neptune having no massive gas envelopes. 

In constructing the simulations presented here, we have appealed to the oligarchic growth 
model as a plausible guide for our choice of initial protoplanet masses and orbits. One can of 
course envision a mimber of variations. For instance, two gas giants may form on initially widely 
separated orbits. The waves they raise in the disk will tend to clear the gas between them, and Type 

11 migration will cause their orbital separation to decrease (Klcy 2000). If ice giant sized bodies 
are trapped in between, they will be prevented from accreting more gas, and will be scattered 
as the stable region between the gas giants shrinks to nothing. Also, if future measurements of 
Jupiter constrain its core to be much smaller, or even absent, this will invalidate the assumption 
that Jupiter grew from an ice giant sized nucleus. Though the scattering of ice giants could still 
take place in principal, one would then need to explain how a small body won the gas accretion 
race against much larger bodies (if the core is small), or how ice giant sized bodies managed to 
form before the birth of the first gas giant from an unnucleated disk instability (if there is no core; 
eg. Boss 2000) . Alternatively, one could search for a way in which accretion could continue to take 
place in close proximity to a gas giant, with scattering delayed until ice giant sized bodies form. 
For such a scenario Type I tidal effects, which cause the eccentricity to decay on an even shorter 
timescalc than the semimajor axis (eg. Papaloizou & Larwood 2000), may actually be helpful. It 
may then be the dispersal of the gas which triggers scattering. Doing away with any reference to a 
specific formation process, the most general statement of our results is: Ice giant sized bodies can be 
scattered from the Jupiter-Saturn region by gas giant sized bodies, to ultimately end up on Uranus- 
and Neptune-like orbits. 
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